##### Extended Data Fig3 - Geographic variation of mutagenic exposures in kidney cancer genomes
library(tidyverse)
library(viridis)
library(reshape2)

#Import Data
#Data in COSMIC fit
data <- read.table(file="ExtendedDataFig3_Attribution.txt",
                   header = TRUE, row.names = 1)

#Transpose so that each sample = 1 column
data <- as.data.frame(t(data))
#Convert Mutation Counts to Relative Contribution
data_relative <- scale(data, center = FALSE, scale = colSums(data))
head(colSums(data_relative))
#Restore Format
data_relative <- as.data.frame(t(data_relative))
data_relative <- tibble::rownames_to_column(data_relative, var="Samples")

#Import Sample Key
sample_key <- read.csv(file="ExtendedDataFig3_Sample_Key.csv")

data_relative=data_relative[order(data_relative$Samples),]
sample_key=sample_key[order(sample_key$Sample),]
data_relative$Cohort <- sample_key$Cohort[sample_key$Sample %in% data_relative$Samples]

#Remove Cohorts with less than 5 cases
data_relative <- filter(data_relative, Cohort != "Cervix-AdenoCA" & Cohort != "Breast-DCIS" & Cohort != "Myeloid-MDS")

#Select Columns to Plot
sbs40_data_relative <- as.data.frame(data_relative[,(c(1,48,47,49,64))])

#Create Proportions
sbs40_proportions <- sbs40_data_relative %>% 
  group_by(Cohort) %>% 
  summarise(
    SBS40a = sum(SBS40a >0)/n(), 
    SBS40b = sum(SBS40b >0)/n(), 
    SBS40c= sum(SBS40c >0)/n(), 
  )

#Summarise (SBS40 Signatures, Relative)
summary <- sbs40_data_relative %>% group_by(Cohort) %>% summarise_if(is.numeric, mean)

#Melt
sbs40_proportions <- melt(sbs40_proportions)
summary_melt <- melt(summary)

#Combine
summary_melt$prop = sbs40_proportions$value

#Plot
#Plot in Landscape with legend at bottom
pdf(file= "ExtendedDataFig3.pdf", height = 2.5, width =7)
ggplot(summary_melt, aes(x = Cohort, y = fct_rev(variable))) +
  geom_rect(ymin = 0, ymax = 4, xmin = 16.5, xmax = 17.5, alpha = 1, color = "black", fill = "grey75", linewidth = 0.25) +
  geom_point(aes(fill = value,  size = prop), color = "black", shape = 21) +
  scale_fill_viridis() +
  theme_bw() +
  theme(axis.title.x = element_blank()) +
  theme(axis.title.y = element_blank()) +
  theme(axis.text.x = element_text(angle = 45, hjust=1, size = 6)) +
  theme(axis.text.y = element_text(size = 6)) +
  labs(size="Percentage of Cases\nwith Signature", fill="Mean Relative\nAttribution") +
  theme(legend.position="bottom") +
  theme(legend.text=element_text(size=6)) +
  theme(legend.title=element_text(size=6))
dev.off()

